4/26/2019

Bayesilaiset Menetelmät

  • Todennäköisyys ja Bayesilainen tilastotiede
    • Bayesin kaava
    • Uskottavuusfunktio
    • Priorijakauma
    • Posteriorijakauma
  • Bayesilaiset menetelmät käytännössä
  • Kalvot löytyvät osoitteesta github.com/hpesonen

Todennäköisyys ja Bayesilainen tilastotiede

  • Bayesiläisen tilastotieteen keskiössä on maailman mallittaminen
  • Bayesilainen tilastotiede mallittaa todennäköisyyden epävarmuutena todellisuudesta
  • 0 ja 1 ovat skaalan kaksi ääripäätä ja vastaavat täyttä varmuutta tapahtumasta
    • Epävarmuus on subjektiivinen käsite
    • Maailman havannointi lisää informaatiota todellisuudesta ja vähentää epävarmuutta
  • Bayesin teoreeman käyttö ei tee kenestäkään Bayesilaista, mutta epävarmuuden ilmaiseminen todennäköisyyden avulla tekee

Todennäköisyys ja Bayesilainen tilastotiede

  • 0 ja 1 ovat skaalan kaksi ääripäätä ja vastaavat täyttä varmuutta tapahtumasta

Todennäköisyys ja Bayesilainen tilastotiede

  • Informaation lisäämisen mekanismi ilmaistaan Bayesin teoreeman avulla
    • Bayesin teoreema on siis työkalu, jota käytetään Bayesilaisessä analyysissä

\[p(\theta \mid y) = \frac{p(y \mid \theta)p(\theta)}{p(y)}\]

Bayesin teoreema

\[p(\theta \mid y) = \frac{p(y \mid \theta)p(\theta)}{p(y)}\]

  • \(\theta\): Mallin parametrit
  • \(y\): Data
  • \(p(\theta \mid y)\): Posteriorijakauma
  • \(p(y \mid \theta)\): Uskottavuusfunktio
  • \(p(\theta)\): Priorijakauma
  • \(p(y)\): Evidenssi/Marginaali

Esimerkki - Kolikonheitto

  • Kolikonheitto on klassinen esimerkki informaation kertymisestä havaintoaineiston kasvaessa
  • Tarkastellaan kolikkoa ja halutaan selvittää onko kolikko reilu
    • Merkitään parametrilla \(p\) todennäköisyyttä saada kolikonheitossa klaava
  • Voidaan ajatella, että kolikonheitto on Bernoullin koe
    • Klaava (\(y=1\)) saadaan todennäköisyydellä \(\theta\), jolloin kruunan todennäköisyys on \(1-\theta\)
    • Uskottavuusfunktio Bernoullin kokeelle on

\[p(y \mid \theta) = \theta^y(1-\theta)^{1-y}\]

Esimerkki - Kolikonheitto

  • Toistettaessa riippumattomasti kolikonheitto-koetta \(N\) kertaa uskottavuusfunktio voidaan ilmaista muodossa

\[p(y_1, \ldots, y_N \mid \theta) = \prod_{i=1}^N p(y_i \mid \theta) = \theta^{\sum_{i} y_i} (1- \theta)^{N-\sum_{i}y_i} \]

Esimerkki - Kolikonheitto

  • Miten mallintaisimme epävarmuuttamme ennen ensimmäistäkään havaintopistettä?
    • Minkälainen priorijakauma vastaisi epävarmuuttamme klaavan esiintymisestä?

Esimerkki - kolikonheitto

  • Uskoaksemme on hyvin epätodennäköistä, että melko varmasti saataisiin aina \(0/1\)
  • Uskoaksemme kolikko on melko todennäköisesti reilu (\(\theta=0.5\)), mutta emme voi poissulkea mahdollisuutta, että se ei ole
  • Mallinnamme prioritietämystä jakaumalla \(p(\theta) \propto \theta^4(1-\theta)^4\) (\(\mathsf{Beta}(5,5)\)-jakauma)

Esimerkki - Kolikonheitto

  • Tarkastellaan miten posteriorijakauma muuttuu heittosarjan koon kasvaessa
  • Merkitään dataa \(y = (n,N-n)\) :
  • Esim. \(y = (10,5)\) vastaa \(15\) heiton sarjaa, joista \(10\) klaavaa ja \(5\) kruunaa
  • Tässä tapauksessa posteriori-jakauma voidaan ilmaista analyyttisesti

\[ \begin{aligned} p(\theta \mid y = (n,N-n)) & \propto p(y=(n,N-n) \mid \theta) p(\theta) \\ & = \theta^{n}(1-\theta)^{N-n}\theta^4(1-\theta)^4 \\ & = \theta^{n+4}(1-\theta)^{N-n+4} \end{aligned} \]

Esimerkki - Kolikonheitto

  • Havaintojen määrän kasvaessa posteriori-jakauma lähestyy uskottavuusfunktiota
    • Priori-todennäköisyyden merkitys vähenee

Esimerkki - Kolikonheitto

  • Tulkitaan posteriorijakauma

Todennäköisyysväli

  • Frekventistisen tilastotieteen luottamusväliä vastaava käsite Bayesilaisessa tilastotieteessa on todennäköisyysväli
  • Määritelmän mukaisesti \(\alpha\)-todennäkäisyysväli (yleisemmin todennäköisyysalue) on mikä tahansa parametrijoukko, joka sisältää \(\alpha\) todennäköisyyttä \[p(\theta \in \Theta \mid x) = \alpha\]
  • Käytännössä kyse on integraalista (tai summasta diskreetin satunnaismuuttujan tapauksessa) \[p(\theta \in \Theta \mid x) = \int_\Theta p(\theta \mid x) ~\mathsf{d} \theta = \alpha\]

Todennäköisyysväli

  • Todennäköisyysvälit eivät ole yksikäsitteisiä
    • Pienin väli/alue, joka sisältää \(\alpha\) todennäköisyyden
    • Väli, jonka alapuolelle jäävä todennäköisyys vastaa yläpuolelle jäävä todennäköisyyttä
    • Odotusarvon ollessa olemassa, väli joka on keskitetty odotusarvoon

Esimerkki - Kolikonheittoesimerkin todennäköisyysväli

  • \(0.95\)-Todennäköisyysväli lyhenee havaintojen määrän kasvaessa

Todennäköisyysväli

  • Ero frekventistiseen luottamusväliin on, että Bayesilaisen tulkinnan mukaisesti parametrin arvo kuuluu \(\alpha\)-todennäköisyysväliin todennäköisyydellä \(\alpha\)

Mallien vertailu

  • Mallien vertailu tapahtuu mallien posterioritodennäköisyyden avulla
  • Jotta mallin \(M_i\) todennäköisyys voidaan laskea, tarvitaan malleille prioritodennäköisyydet \(p(M_i)\) kuten mille tahansa parametrille

\[\sum_{i=1}^m p(M_i) = 1\]

  • Mallien posterioritodennäköisyydet ovat

\[\sum_{i=1}p(M_i \mid y) = 1\]

Mallien vertailu

  • Oletetaan, että mallin \(M_i\) määrittää sen parametrit \(\theta_i\)
  • Tällöin mallin posterioritodennäköisyys lasketaan \[ p(M_i \mid y) = \frac{p(y \mid M_i) p(M_i)}{p(y)} \]
  • Nyt mallin uskottavuus, evidenssi, lasketaan \[p(y \mid M_i) = \int p(y \mid \theta_i, M_i) ~ p(\theta_i \mid M_i) ~ \mathsf{d} \theta_i\]

Mallien vertailu

  • Nykyisin käytetään usein mallien vertailuun posteriorijakauman laskemisen jälkeen Leave-one-out ristiinvalidointia
  • Laskettava suure on \[\text{elpd}_\text{loo} = \sum_{i=1}^N \log p(y_i \mid y_{-i}) = \sum_{i=1}^N \log \int p(y_i \mid \theta)p(\theta \mid y_{-i}) \mathsf{d}\theta\]

Prediktiiviset jakaumat

  • Suureen \(\text{elpd}_\text{loo}\) laskemissa käytetään hyväksi prediktiivisiä jakaumia
  • Prediktiivisiä jakaumia voidaan yleisemminkin hyödyntää havaintomallien sekä priorijakaumien järkevyyden tarkasteluun
  • Nämä voidaan jakaa prioriprediktiivisiin sekä posterioriprediktiivisiin jakaumiin

Prioriprediktiivinen jakauma

  • Kun havaintomalli sekä tämän parametrien priorijakauma on mallitettu, minkälaista havaintoaineistoa malli ennustaa

\[p(y) = \int p(y \mid \theta) p(\theta) \mathsf{d} \theta\]

  • Usein prioritieto liittyy nimenomaisesti tietämykseen mallin tuottamista realistisista arvoista

  • Prioriprediktiivinen malli selventää parametrien priorijakauman vaikutusta havantoihin
    • Tuottavatko jotkut parametrikombinaatiot mahdottomia/hyvin epätodennäköisiä havaintoja?

Posteoriprediktiivinen jakauma

  • Kun havaintomalli sekä tämän parametrien priorijakauma on mallitettu ja posteriorijakauma on laskettu havaintoaineistoa käyttäen, minkälaista uutta havaintoaineistoa \(y^\prime\) malli ennustaa

\[ \begin{aligned} p(y^\prime \mid y) & = \int p(y^\prime \mid \theta, y) p(\theta \mid y) \mathsf{d} \theta \\ & = \int p(y^\prime \mid \theta) p(\theta \mid y) \mathsf{d} \theta \end{aligned} \]

  • Posterioriprediktiivistä mallinnusta voidaan käyttää tarkastelemaan opittujen mallien järkevyyttä

Bayesilainen analyysi käytännössä

  • Posteriori-jakauman laskeminen analyyttisessä muodossa onnistuu ainoastaan erikoistapauksissa
  • Käytännössä laskeminen tapahtuu muodostamalla otosparvi posteriorijakaumasta
    • Importance sampling
    • Sequential Monte Carlo
    • Markov Chain Monte Carlo
    • Hamiltonian Monte Carlo

Otosparvi posteriorijakauman estimaattina

  • Otosparven avulla voidaan likimääräisesti laskea \(\mathsf{E}(g( \theta)\mid y)\) mille tahansa funktiolle \(g(\cdot)\)

Stan

  • Nykyisin otoksen muodostamiseen on yksinkertaista laajaltikäytettyjen ohjelmistojen ansiosta
  • Suosituin ohjelmistoista on Stan
    • Stan
      • PyStan
      • RStan
  • rstanarm : Bayesian applied regression modeling
    • R:stä tuttu käyttöliittymä Bayesilaisesta analyysista kiinnostuneille

rstanarm-esimerkki

  • Data : kidiq
  • Muuttujat
    • Lapsen ÄO
    • Äidillä high school-tutkinto indikaattorimuuttuja
    • Äidin ÄO
    • Äidin ikä

rstanarm-esimerkki

rstanarm-esimerkki

Esimerkki

Esimerkki

Esimerkki

Esimerkki 2

Esimerkki 3

Esimerkki 4

Esimerkki

Esimerkki

Esimerkki

## 
## Model comparison: 
## (ordered by highest ELPD)
## 
##       elpd_diff se_diff
## post4   0.0       0.0  
## post3  -3.5       2.5  
## post2  -6.0       3.9  
## post1 -42.2       8.7

Esimerkki

## 
## Model comparison: 
## (negative 'elpd_diff' favors 1st model, positive favors 2nd) 
## 
## elpd_diff        se 
##      42.2       8.7

Esimerkki

## 
## Model comparison: 
## (negative 'elpd_diff' favors 1st model, positive favors 2nd) 
## 
## elpd_diff        se 
##       3.5       2.5
## 
## Model comparison: 
## (negative 'elpd_diff' favors 1st model, positive favors 2nd) 
## 
## elpd_diff        se 
##       6.0       3.9

Esimerkki

Esimerkki

Esimerkki

Esimerkki

## 
## Call:
## glm(formula = lot1 ~ log(u), family = Gamma, data = clotting)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -0.04008  -0.03756  -0.02637   0.02905   0.08641  
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0165544  0.0009275  -17.85 4.28e-07 ***
## log(u)       0.0153431  0.0004150   36.98 2.75e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 0.002446059)
## 
##     Null deviance: 3.51283  on 8  degrees of freedom
## Residual deviance: 0.01673  on 7  degrees of freedom
## AIC: 37.99
## 
## Number of Fisher Scoring iterations: 3
## 
## Call:
## glm(formula = lot2 ~ log(u), family = Gamma, data = clotting)
## 
## Deviance Residuals: 
##      Min        1Q    Median        3Q       Max  
## -0.05574  -0.02925   0.01030   0.01714   0.06371  
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.0239085  0.0013265  -18.02 4.00e-07 ***
## log(u)       0.0235992  0.0005768   40.91 1.36e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 0.001813354)
## 
##     Null deviance: 3.118557  on 8  degrees of freedom
## Residual deviance: 0.012672  on 7  degrees of freedom
## AIC: 27.032
## 
## Number of Fisher Scoring iterations: 3

Esimerkki

Esimerkki

## stan_glm
##  family:       Gamma [inverse]
##  formula:      clot_time ~ log_plasma * lot_id
##  observations: 18
##  predictors:   4
## ------
##                    Median MAD_SD
## (Intercept)        -0.016  0.007
## log_plasma          0.015  0.003
## lot_id2            -0.007  0.013
## log_plasma:lot_id2  0.008  0.006
## 
## Auxiliary parameter(s):
##       Median MAD_SD
## shape 7.942  2.577 
## 
## Sample avg. posterior predictive distribution of y:
##          Median MAD_SD
## mean_PPD 31.944  4.859
## 
## ------
## * For help interpreting the printed output see ?print.stanreg
## * For info on the priors used see ?prior_summary.stanreg

Esimerkki

Slide with R Output

Esimerkki

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Esimerkki

Esimerkki

Esimerkki

## (Intercept)     dist100     arsenic 
##       0.000      -0.894       0.461

Esimerkki

Esimerkki

## Loading required package: gridExtra
## 
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
## 
##     combine
## Warning: Removed 39 rows containing missing values (geom_point).